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Abstract 

We present a new algorithm for reducing an arbitrary unitary matrix into a sequence 
of elementary operations (operations such as controlled- nots and qubit rotations). 
Such a sequence of operations can be used to manipulate an array of quantum bits 
(i.e., a quantum computer). We report on a C++ program called "Qubiter" that 
implements our algorithm. Qubiter source code is publicly available. 
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1. Introduction 



1(a) Previous Work 

In classical computation and digital electronics, one deals with sequences of 
elementary operations (operations such as AND, OR and NOT). These sequences 
are used to manipulate an array of classical bits. The operations are elementary in 
the sense that they act on only a few bits (1 or 2) at a time. Henceforth, we will 
sometimes refer to sequences as products and to operations as operators, instructions, 
steps or gates. Furthermore, we will abbreviate the phrase "sequence of elementary 
operations" by "SEO". In quantum computation 0, one also deals with SEOs (with 
operations such as controlled- nots and qubit rotations), but for manipulating quan- 
tum bits (qubits) instead of classical bits. SEOs are often represented graphically by 
a qubit circuit. 

In quantum computation, one often knows the unitary operator U that de- 
scribes the evolution of an array of qubits. One must then find a way to reduce U 
into a SEO. In this paper, we present a new algorithm for accomplishing this task. 
We also report on a C++ program called "Qubiter" that implements our algorithm. 
Qubiter source code is publicly available at www.ar-tiste.com/qubiter.html . We call 
Qubiter a "quantum compiler" because, like a classical compiler, it produces a SEO 
for manipulating bits. 

Our algorithm is general in the sense that it can be applied to any unitary 
operator U. Certain U are known to be polynomial in Nb] i-e, they can be expressed 
as a SEO whose number of gates varies as a polynomial in Nb as Nb varies, where 
Nb is the number of bits. For example, a Discrete Fourier Transform (DFT) can be 
expressed as 0rder(iV5 2 ) steps 0. Our algorithm expresses a DFT as Order (4^) 
steps. Hence, our algorithm is not "polynomially efficient"; i.e., it does not give 
a polynomial number of steps for all U for which this is possible. However, we 
believe that it is possible to introduce optimizations into the algorithm so as to 
make it polynomially efficient. Future papers will report our progress in finding such 
optimizations. 

Previous workers |3j have described another algorithm for reducing a unitary 
operator into a SEO. Like ours, their algorithm is general, and it is not expected to 
be polynomial efficient without optimizations. Our algorithm is significantly different 
from theirs. Theirs is based on a mathematical technique described in Refs.@j-[|J, 
whereas ours is based on a mathematical technique called CS Decomposition|J-[^|, 
to be described later. 

Quantum Bayesian (QB) Nets[B]-[]IOj are a method of modelling quantum sys- 



tems graphically in terms of network diagrams. In a companion paper pi]] , we show 
how to apply the results of this paper to QB nets. 
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1(b) CS Decomposition 



As mentioned earlier, our algorithm utilizes a mathematical technique called 
CS Decomposition^]- j5|. In this name, the letters C and S stand for "cosine" and 
"sine" . Next we will state the special case of the CS Decomposition Theorem that 
arises in our algorithm. 

Suppose that U is an iV x N unitary matrix, where N is an even number. 
Then the CS Decomposition Theorem states that one can always express U in the 
form 



U 



where L , Li, R , R\ are y x y un itary matrices and 
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D 00 = Dn = diag(Ci, C 2 , . . . , Cn) 
diag(Si, S 2 , ■ ■ ■ , Sn) , 
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For i e {1,2, 



, y}, Ci and Si are real numbers that satisfy 



Cf + Sf 



(1.2a) 
(1.2b) 
(1.2c) 
(1.2d) 

(1.2e) 



Henceforth, we will use the term D matrix to refer to any matrix that satisfies 
Eqs.(1.2). If one partitions U into four blocks Uij of size y x y , then 



(1.3) 

for i, j G {0, 1}. Thus, gives the singular values of U^. 

More general versions of the CS Decomposition Theorem allow for the possi- 
bility that we partition U into 4 blocks of unequal size. 

Note that if U were a general (not necessarily unitary) matrix, then the four 
blocks Uij would be unrelated. Then to find the singular values of the four blocks Uij 
would require eight unitary matrices (two for each block), instead of the four Li, Rj. 
This double use of the Lj, Rj is a key property of the CS decomposition. 
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1(c) Bird's Eye View of Algorithm 

Our algorithm is described in detail in subsequent sections. Here we will only 
give a bird's eye view of it. 




L(00, UL) R(00, UL) 

L(01, UL) D(0, UL) R( oi, UL) 
L(10,UL) D(1,UL) R(10,UL) 
L(ll, UL) R(ll, UL) 




L(00, UR) R( oo, UR) 

L(01,UR) D(0,UR) R(01,UR) 
L(10, UR) D(1,UR) R(10,UR) 
L(ll, UR) R(ll, UR) 




Fig.l A binary tree. Each node (3 has a single parent. If the parent is to (3's right 
(ditto, left), then (3 contains the names of the matrices produced by applying the 
CS Decomposition Theorem to the L matrices (ditto, R matrices) of /3's parent. 

Consider Fig.l. We start with a unitary matrix U. Without loss of generality, 
we can assume that the dimension of U is 2 Nb for some Nb > 1. (If initially U's 
dimension is not a power of 2, we replace it by a direct sum U @ I T whose dimension 
is a power of two.) We apply the CS Decomposition method to U. This yields a 
matrix D(0, U) of singular values, two unitary matrices L(0, U) and L(l, U) on the 
left and two unitary matrices R(0, U) and R(l, U) on the right. Then we apply the CS 
Decomposition method to each of the 4 matrices L(0, U), L(l, U), R(0, U) and R(l, U) 
that were produced in the previous step. Then we apply the CS Decomposition 
method to each of the 16 R and L matrices that were produced in the previous step. 
And so on. The lowest row of the pyramid in Fig.l has L's and i?'s that are lxl 
dimensional, i.e., just complex numbers. 
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Call a central matrix either (1) a single D matrix, or (2) a direct sum D\ © 
D 2 © • • • © D r of D matrices, or (3) a diagonal unitary matrix. From Fig.l it is clear 
that the initial matrix U can be expressed as a product of central matrices, with each 
node of the tree providing one of the central matrices in the product. Later on we 
will present techniques for decomposing any central matrix into a SEO. 



2. Preliminaries 

In this section, we introduce some notation and some general mathematical 
concepts that will be used in subsequent sections. 

2(a) General Notation 

We define Z a ^ = {a, a + 1, . . . , b} for any integers a and b. S(x, y) equals one 
if x = y and zero otherwise. 

We will use the symbol Nb for the number (> 1) of bits and N$ = 2 Nb 
for the number of states with Nb bits. Let Bool = {0, 1}. We will use lower case 
Latin letters a,b,c . . . G Bool to represent bit values and lower case Greek letters 
a, (3, 7, . . . G Z 0i jv b _i to represent bit positions. A vector such as a = a N[3 _i . . . a 2 aia 
will represent a string of bit values, a M being the value of the /x'th bit for \x G Z 0i jv s -i- 
A bit string a has a decimal representation d(a) = X^Jo" 1 2 fl a fl . For j3 G Zq^-i, we 
will use u(P) to denote the /3'th standard unit vector, i.e, the vector with bit value of 
1 at bit position (3 and bit value of zero at all other bit positions. 

I r will represent the r dimensional unit matrix. Suppose (3 G Z 0i Ar B _i and M 
is any 2x2 matrix. We define M(j3) by 

M(J3) = h © • • • © h ® M © h © • • • © h , (2.1) 

where the matrix M on the right side is located at bit position f3 in the tensor product 
of Nb 2x2 matrices. The numbers that label bit positions in the tensor product 
increase from right to left (■*—), and the rightmost bit is taken to be at position 0. 

For any two square matrices A and B of the same dimension, we define the 
product by A © B = ABA\ where A^ is the Hermitian conjugate of A. 

a = (o~ x , o~ y ,o~ z ) will represent the vector of Pauli matrices, where 





( 






< -i \ 




o- x = 




o) 


, (Jy = 


v i J 






2(b) Sylvester-Hadamard Matrices 

The Sylvester-Hadamard matrices |E| H r are defined by: 
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(2.3a) 



H r+ i = Hi (g) H r 



(2.3b) 



(2.3c) 



for any integer r > 1. In Eqs.(2.3), we have labelled the rows and columns with 
binary numbers in increasing order, and ® indicates a tensor product of matrices. 
From Eqs.(2.3), one can show that the entry of H r at row a and column b is given by 



\a-b 



where a ■ b = X^ =0 a^b^. It is easy to check that 



Hi 



Hi 



H r , 



2 r L 



(2.4) 

(2.5) 
(2.6) 



In other words, H r is a symmetric matrix, and the inverse of H r equals H r divided 
by however many rows it has. 

2(c) Permutations 

Subsequent sections will use the following very basic facts about permutations. 



For more details, see, for example, Ref. . 

A permutation is a 1-1 onto map from a finite set X onto itself. The set of 
permutations on set X is a group if group multiplication is taken to be function 
composition. S n , the symmetric group in n letters, is defined as the group of all 
permutations on any set X with n elements. If X = Zi_ n , then a permutation G 
which maps i G X to cij G X (where i ^ j implies a, 7^ aj) can be represented by a 
matrix with entries 



(G)j,i = S(aj,i) 



(2.7) 



for all i,j G X. Note that all entries in any given row or column equal zero except 
for one entry which equals one. An alternative notation for G is 
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G= 1 2 3 n ) . (2.8) 

\ a ± a 2 a 3 ■ ■ ■ a n J 

The product of two symbols of the type shown in Eq.(2.8) is defined by function 
composition. For example, 

fa 1 a 2 a 3 \fl2 3\fl2 3\ 
\ h b 2 b 3 )\ ai a 2 a 3 ) \ bi b 2 b 3 J 

Note how we have applied the permutations on the left side from right to left (<—). 
(Careful: Some authors apply them in the opposite direction (—>)). A cycle is a 
special type of permutation. If G e S n maps oi — > a 2 , a 2 — > a 3 , . . . , a r _i — > a r , 
a r — > ai, where « 7^ j implies a, 7^ and r < n, then we call G a cycle. G may be 
represented as in Eqs.(2.7) and (2.8). Another way to represent it is by 

G = (ai,a 2 ,a 3 , ...,a r ). (2.10) 

(Careful: some people write (a r , . . . , a 3 , a 2 , ai) instead.) We say that the cycle of 
Eq.(2.10) has length r. Cycles of length 1 are just the identity map. A cycle of length 
2 is called a transposition. The product of two cycles need not be another cycle. For 
example, 

(2, 1,5)(1,4,5,6) = (I 1 3 2 6 5) (2 ' n) 

cannot be expressed as a single cycle. Any permutation can be written as a product 
of cycles. For example, 

(I 1 3 2 6 S)=(MXM,2). (2.12) 

The cycles on the right side of Eq.(2.12) are disjoint; i.e., they have no elements 
in common. Disjoint cycles commute. Any cycle can be expressed as a product of 
transpositions (assuming a group with > 2 elements), by using identities such as: 

(ai, a 2 , . . . , a n ) = (ai, a 2 )(a 2 , a 3 ) ■ ■ ■ (a„_i, a n ) , (2.13) 

(ai, a 2 , . . . , a„) = (ai, a„) • • • (a u a 3 )(ai, a 2 ) . (2.14) 
Another useful identity is 

(a,6) = (a,p)(p,6)(a,p). (2.15) 
This last identity can be applied repeatedly. For example, applied twice, it gives 
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(a, b) = (a,pi)(pi, 6)(a,pi) = (a,p 1 )(p 1 ,p 2 )(p2, b){p 1 ,p 2 )(a,p 1 ) . (2.16) 



Since any permutation equals a product of cycles, and each of those cycles can be 
expressed as a product of transpositions, all permutations can be expressed as a prod- 
uct of transpositions (assuming a group with > 2 elements). The decomposition of 
a permutation into transpositions is not unique. However, the number of transposi- 
tions whose product equals a given permutation is always either even or odd. An even 
(ditto, odd) permutation is defined as one which equals an even (ditto, odd) number 
of transpositions. 



2(d) Projection Operators 

Consider a single qubit first. 

The qubit's basis states |0) and |1) will be represented by 



10} 



1 






1 



The number operator n of the qubit is defined by 



n 




1 



1 — 0", 



Note that 



(2.17) 



(2.18) 



n\0) = , n\l) = |1) 
We will often use n as shorthand for 



n = 1 — n 



1 




1 + 0-* 



Define Pq and Pi by 



P n = n 



1 




P 1 = n 




1 



Po and Pi are orthogonal projectors and they add to one: 

P a P b = S(a, b)P b for a, b e Bool 



(2.19) 



(2.20) 



(2.21} 



(2.22a) 



P + Pi = h ■ (2.22b) 
Now consider iV^ bits instead of just one. 

For (3 e Z , Nb -i, we define P (/3), Pi{fl), n{0) andn(/3) by means of Eq. (2.1). |4] 
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For ae BooI Nb , let 

Pa = Pa NB ., ®---®Pa 2 ®P ai ® P ao • (2.23) 

Note that 

Pa = h ® h ® ■ ■ ■ ® h = I 2 »b • (2.24) 

a<=Bool N B 

For example, with 2 bits we have 

P O0 = P ®P = diag(l,0,0,0), (2.25a) 

P 01 = P ®P 1 = diag(0, 1, 0, 0) , (2.25b) 

P 10 = p 1 ®P = diag(0, 0,1,0) , (2.25c) 

P u = Pi (8) Pi = diag(0, 0, 0, 1) , (2.25d) 

E ^, 6 = / 4 . (2.26) 

For r > 1, suppose P 1 ,P 2 ,...P r are orthogonal projection operators (i.e., 
pPj = 5(i,j)Pj ), and cti , «2 . . . a r are complex numbers. Then it is easy to show by 
Taylor expansion that 

r r r 

exp(]T otiPi) = E exp( ai )P + (1 - E ^) • (2-27) 

i=l i=l i=l 

In other words, one can "pull out" the sum over Pj's from the argument of the 
exponential, but only if one adds a compensating term 1 — ^ Pi so that both sides 
of the equation agree when all the a^s are zero. 



3. State Permutations that Act on Two Bits 

The goal of this paper is to reduce unitary matrices into qubit rotations and 
controlled-nots (c-nots). A qubit rotation (i.e., exp[i9 ■ <?(/?)] for (5 G Zq^ b ^i and 
some real 3-dimensional vector 9) acts on a single qubit at a time. This section will 
discuss gates such as c-nots that are state permutations that act on two bits at a 
time. 
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3(a) N B = 2 

Consider first the case when there are only 2 bits. Then there are four possible 
states-00, 01, 10, 11. With these 4 states, one can build 6 distinct transpositions: 



(00,01) 
(00, 10) 



(01,11) 



(10,11) 



Pi Po 
Po Pi 



(00,11) 



(01,10) 



Po Pi 
Pi Po 



a x 



P ®a x + P l ®I 2 = a x (0) n(1) 

h®Pi + a x ®P a = a x {lf^ 
1 



1 
1 



h ® Po + <T X ® Pi = ^(1) 



n(0) 



P ®I 2 + Pi®a x = cr x (0) 



n(l) 



(3.1a) 
(3.1b) 

(3.1c) 

(3.1d) 

(3.1e) 
(3-lf) 



where matrix entries left blank should be interpreted as zero. The rows and columns of 
the above matrices are labelled by binary numbers in increasing order (as in Eq.(2.3b) 
for H 2 ). Note that the 4 transpositions Eqs.(3.1 a,b,e,f) change only one bit value, 
whereas the other 2 transpositions Eqs.(3.1 c,d) change both bit values. We will call 
(00,11) the Twin-to-twin- er and (01,10) the Exchanger. 



Exchanger [jig] has four possible representations as a product of c-nots: 



(oi, io) = (oi,oo)(oo,io)(oi,oo) = a x (or(V(i) n(c V(or (1) , 



(3.2a) 



(01, 10) = (10, 11)(11, 01)(10, 11) = a x (0)"(V(ir ( V(0r (1) 



(3.2b) 



(01,10) = (io,oo)(oo,oi)(io,oo) = ^(ir(V(or ( V(i) n(0) , 



(3.2c) 



(01, 10) = (01, 11)(11, 10)(01, 11) = a x (ir^a x (0r^a x (ir^ . 



(3.2d) 
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Note that one can go from Eq.(3.2a) to (3.2b) by exchanging n and n; from Eq.(3.2a) 
to (3.2c) by exchanging bit positions and 1; from Eq.(3.2a) to (3. 2d) by doing both, 
exchanging n and n and exchanging bit positions and 1. We will often represent 
Exchanger by E(0, 1). It is easy to show that 

E T (0, 1) = £(0,1) = £ _1 (0, 1) , (3.3a) 
£(0,1) = E(1,0) , (3.3b) 

£ 2 (0,1) = 1. (3.3c) 

Furthermore, if X and Y are any two arbitrary 2x2 matrices, then, by using the 
matrix representation Eq.(3.1d) of Exchanger, one can show that 

E(l^)Q(X®Y) = Y®X . (3.4) 

Thus, Exchanger exchanges the position of matrices X and Y in the tensor product. 
Twin-to-twin-er also has 4 possible representations as a product of c-nots. One 

is 

(00,11) = (00,01)(01,ll)(00,01) = ( r ;c (0f( 1 )a :c (ir( )(T :c (0f« . (3.5) 

As with Exchanger, the other 3 representations are obtained by replacing (1) n and 
n, (2) bit positions and 1, (3) both. 
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a (0)n(l)= k — ± 

A. 



a x (l) 



n(0) = 



X <? 



a x (l) n (°) = >^ 



t 



■x 



a x (0) n (D 



Fig. 2 Circuit symbols for the 4 different types of c- riots. 







(01,10) = X 6 = X 



r 



y 



x x 



A 



t 



x 



<> — x 
x^ 



Fig. 3 Four equivalent circuit diagrams for Exchanger. 
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(00,11) = ± 



y 



y 



= X 



y 



A 



o - 



y 



x 



X = 



o- 

X 



Fig. 4 Four equivalent circuit diagrams for Twin-to- twin-er. 

Figures 2, 3 and 4 give a diagrammatic representation of the 6 possible trans- 
positions of states for Nb = 2. 

3(b) Any N B > 2 

Suppose ai, &i, a 2 , b 2 G Boo/ and a, (3 & Z n b -i- We define 



(a 1 fei,a 2 6 2 ] 



o,/3 



n 



(AaiA'ftiA", Aa 2 A / 6 2 A") 



(3.6) 



(A,A' ,A")&BooV 



where on the right side, a 1; a 2 are located at bit position a, and &i, 5 2 ar e located at 
bit position (3. (Note that the transpositions on the right side of Eq.(3.6) are disjoint 
so they commute.) For example, for Nb = 3, 





4 

y 



t 



^(O)^ = (10,ll)i,o= II 010, all) = (010, 011)(110, 111) . (3.7) 

adBool 

Clearly, any permutation of states with Nb bits that permutes only 2 bits (i.e., Ex- 
changer, Twin-to-twin-er, and all c-nots) can be represented by (ai&i, a2& 2 )a,/3- 
For a, (3 G Z 0j n b -i, let E(a,/3) represent Exchanger: 

E(a,/3) = (01,10)^. (3.8) 

As in the Nb — 2 case, E(a, (3) can be expressed as a product of c-nots in 4 different 
ways. One way is 



(3.9) 



The other 3 ways are obtained by exchanging (1) n and n, (2) bit positions a and /3, 
(3) both. Again as in the Nb = 2 case, 



E T {a,P) = E(a,(3) = E^(a,(3) 



(3.10a) 



E(a,{3) = E(f3,a) 



(3.10b) 
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E 2 (a,P) = l. 



(3.10c) 



Furthermore, if X and Y are two arbitrary 2x2 matrices and a,/3 G Zo,Ar s -i such 
that a 7^ (3, then 



S(a,/3) [X(a)F 



X(/3)Y(a) . 



(3.11) 



Equation (3.11) is an extremely useful result. It says that E(a, (3) is a transpo- 
sition of bit positions. (Careful: this is not the same as a transposition of bit states.) 
Furthermore, the E(a,(3) generate the group of Nb\ permutations of bit positions. 
(Careful: this is not the same as the group of (2^)! permutations of states with Nb 
bits.) 



E(a, P) 



= 1 



(3 a 



t T 

Fig. 5 Circuit symbol for Exchanger. 
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Fig. 6 Circuit diagram for Eq.(3.12). 
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2 1 



Fig. 7 Circuit diagram for Eq.(3.14). 
An example of how one can use Eq.(3.11) is 



E(2, 0) [a x (0r^a x (2) n ^a x (3r^} = a x (2) n W a x {0) n ^ a x (3) n ^ . (3.12) 

Figure 5 gives a convenient way of representing E(a, (3) diagrammatically. Using this 
symbol, the example of Eq.(3.12) can be represented by Fig. 6. 

Of course, identities that are true for a general transposition are also true for 
E(a,p). For example, 



(2,0) = (2,1)(1,0)(2,1) 



Therefore, 



£(2,0) = £(2,1)£(1,0)£(2,1). 
Figure 7 is a diagrammatic representation of Eq.(3.14). 



(3.13) 
(3.14) 



4. Decomposing Central Matrix into SEO 

In Section 1(d), we gave only a partial description of our algorithm. In this 
section, we complete that description by showing how to decompose each of the 3 
possible kinds of central matrices into a SEO. 



4(a) When Central Matrix is a Single D Matrix 

D matrices are defined by Eqs.(1.2). They can be expressed in terms of pro- 
jection operators as follows: 



D 



E 



exp(i(p a a. 



id 



Pn 



(4.1) 



where the 0^ are real numbers. Note that in Eq.(4.1), a has Nb — 1 components 
instead of the full N B . Using the identity Eq.(2.27), one gets 
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D = eX P [l J2 <t>aVy®Pa\ ■ (4.2) 
aeBoo^B- 1 



Now define new angles #g by 



a ~ 

beBoo^B- 1 



E (~lf% ■ (4-3) 



Suppose (f) (ditto, 9) is a column vector whose components are the numbers (fig (ditto, 
0s) arranged in order of increasing a. Then Eq.(4.3) is equivalent to 

$=H Nb -J. (4.4) 

This is easily inverted to 

= ^ztHn b -J- (4.5) 
Let A z for b G BooI Nb ~ 1 be defined by 

^=exp(<^® £ (-lf b Ps) ■ (4.6) 

\ aeBoo^B- 1 / 

Then D can be written as 

D= J] A b- (4-7) 

beBoo^B- 1 

Note that the operators on the right side commute so the order in which they 
are multiplied is irrelevant. Next we establish 2 useful identities: If (5 G Z 0) jv s -i and 
«(/?) G BooI Nb is the /3'th standard unit vector, then 



E^Boo^(-l)^ (/3) ^ 



= I 2 <8> ■ ■ ■ ® h € 

= P (/?)-A(/5) 
If /3, a G Z 0i 7v B -i and a 7^ /3, then 



(4.8) 



a*(/?) n(Q) 0<7„(/?) 

= [^(/?)Pi(a) + P o (a)]0a- y (/?) 

= a„(/3)[-P 1 (a) + P (a)] 
= <r v (P)<r z (oi) 



(4.9) 



Now we are ready to express Ag- in terms of elementary operators. For any b G 
Bool 1 ^ 13 ^ 1 , we can write 
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where 



3=0 



(4.10) 



Wjj - 2 > > • • • > ft > A, > . 



(4.11) 



In other words, b has bit value of 1 at bit positions (3j. At all other bit positions, b 
has bit value of 0. r is the number of bits in b whose value is 1. When b = 0, r = 0. 
Applying Eqs.(4.8) and (4.9), one gets 



[a x (N B - • • • a x (N B - l) n ^a x (N B - l) n ^\ a y (N B - 1) 

= (7 W (jVfl-i)n5i5 (7,(^0 

= IlJ=o (EaC-l)* 3 ^^ ® ^a) 

= Ea(-l)"X®P S 

Thus, 



A; = M^b - l) n( ^- l} • • • <r x (N B - l) n ^a x (N B - 1)^] exp[i^ y (JV fl - 1)] . 



A- h = [a x {P r -2) n ^ ■ ■ ■ a x {Pi) nm a x {f3 ) n(M } exp^^A))] . (4.14) 



In conclusion, we have shown how to decompose a D matrix into a SEO. For 
example, suppose N B = 3. Then 




D = E exp(i0 ab cr y ) P a P, . 



(4.15) 



Define # by 




(4.16) 



By Eqs.(4.7) and (4.13), 



D = A 00 A 01 A 10 A n , 



(4.17) 



where 



A 00 = exp(i9 00 a y ) J 2 h , 



(4.18a) 
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A m = a x (2) n ^ © [exp(ie 01 <T y ) J 2 / 2 ] , 



(4.18b) 



Aio = ^(2)^) [exp(i0io<7 y ) / 2 ® / 2 ] , (4.18c) 
A u = [a x (2) n ^a x (2) n ^\ [exp(ie u a y ) I 2 J 2 ] . (4.18d) 

4(b)When Central Matrix is a Direct Sum of D Matrices 

Consider first the case N B = 3. Let R{<f>) = exp(ia y 4>). Previously we used 
the fact that any D matrix D can be expressed as 

D= E R{€ b )®Pa®P b . (4.19) 

a,b£Bool 

But what if R were located at bit positions or 1 instead of 2? By expressing both 
sides of the following equations as 8 x 8 matrices, one can show that 

D ®D 1 = E Pa®R(<P' ab )®Pb, (4.20) 
a,b£Bool 

D 00 © D 01 © D w © D u = E P a®Pb® R((f>ab) , (4.21) 

a,b£Bool 

where the Dj and are D matrices. One can apply a string of Exchangers to move 
R in Eqs.(4.20) and (4.21) to any bit position. Thus, 

D @D 1 = E{1, 2)0 I E R{4>'ab)®Pa®Pb\ , (4.22) 

\a,beBooi / 

A)0 © A)l © #10 © £>11 = £(O,1)£(1,2)0 I E R^ab)®Pa®Pb) ■ (4.23) 

(Careful: £(0,2) ^ £(0, 1)£(1, 2). £(0,2) will change Y,a,bR{4>ab) ® P a ® Pb to 
Sa,b Pb® P a ® R{4>ab)-, which is not the same as the left side of Eq.(4.21) ). 
For general N B > 1, if k e Z 0) jv b -i and 




1 if fc = or iV B = 1 

E(N B - k - 1, N B - k) ■ ■ ■ E(N B -3,N B - 2)E(N B -2,N B - 1) otherwise 



(4.24) 

then a direct sum of 2 k D matrices can be expressed as 
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E I J2 R(.<f>s) ®Ps) ■ (4.25) 

It follows that if we want to decompose a direct sum of D matrices into a SEO, 
we can do so in 2 steps: (1) decompose into a SEO the D matrix that one obtains 
by moving the qubit rotation to bit position Nb — 1, (2) Replace each bit name in 
the decomposition by its "alias" . By alias we mean the new name assigned by the bit 
permutation E defined by Eq.(4.24). 

4(c) When Central Matrix is a Diagonal Unitary Matrix 

Any diagonal unitary matrix C can be expressed as 

C= Y, exp(# s )P ff , (4.26) 
deBooi N B 

where the 4>s are real numbers. Using the identity Eq.(2.27) yields 

C = exp(* ]T faPA . (4.27) 

\ a£Bool N B J 

Now define new angles 9^ by 

fa= £ (~lf%. (4.28) 



beBooi"B 



In terms of vectors, 



and 



2 N B 

Let A? for b E BooI Nb be defined by 



H Nb , (4.29) 



6 = -Lh n J. (4.30) 



*i = «pK E (-lf 6 ^| • (4-31) 

SeBool N B 



Then C can be written as 

C= n A". ( 4 - 32 ) 

beBool N B 

where the operators commute. For any b G BooI Nb , we can write 
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3=0 



(4.33) 



where 



#5 - 1 > P r -! > ■ ■ ■ > A > A) > ■ 



(4.34) 



(Careful: Compare this with Eq.(4.11). Now b E BooI Nb instead of BooI Nb ~ 1 and 
(3 r -i can be as large as Nb — 1 instead of Nb — 2.) One can show using the techniques 
of Section 4(a) that 

exp[z# ] if t = 
exp[i9pr z (Po)] if r = 1 

~<7x(Po) n(l3r - l) ■ ■ ■ *x(l3o) n{fh) <Tx(Po) nm ] © exp[i^(A))] if r > 2 

(4-35) 

As in Section 4(a), there are other ways of decomposing into a SEO. 

In conclusion, we have shown how to decompose a diagonal unitary matrix 
into a SEO. For example, suppose N B = 2. Then 



C = diagie** 00 , e** 01 , e** 10 , e^ 11 ) . 



Define by 



By Eqs.(4.32) and (4.35), 



where 



= ~#20 • 



C — ^4-oo^-oi ^4io^4ii i 



(4.36) 



(4.37) 



(4.38) 



A 00 = exp(i#oo) 



(4.39a) 



An = h® exp(i9 ia z ) 



Aio = exp(i6»i a 2 ) ® I 2 , 



A n = a x ( y 0) n{1) [I 2 ® exp(^ n a 2 



(4.39b) 
(4.39c) 
(4.39d) 
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5. Qubiter 



At present, Qubiter is a very rudimentary program. We hope that its fans will 
extend and enhance it in the future. Qubiterl.O is written in pure C++, and has 
no graphical user interface. In its "compiling" mode, Qubiter takes as input a file 
with the entries of a unitary matrix and returns as output a file with a SEO. In its 
"decompiling" mode, it does the opposite: it takes a SEO file and returns the entries 
of a matrix. The lines in a SEO file are of 4 types: 

(a) PHAS ang 

where ang is a real number. This signifies a phase factor exp (i (ang) ^). 

(b) CNOT a char (3 

where a, (3 G Z ^ B _i and char G {T, F}. T and F stand for true and false. 
If char is the character T, this signifies a x ((3) n ( a \ Read it as "c-not: if a is 
true, then flip (3. n If char is the character F, this signifies cr x (/3) n ( a \ Read it as 
"c-not: if a is false, then flip f3." 

(c) ROTY a ang 

where a G Zq^ b -i and ang is a real number. This signifies the rotation 
of qubit a about the Y axis by an angle ang in degrees. In other words, 
exp(ia y (a)ang^). 

(d) ROTZ a ang 

This is the same as (c) except that the rotation is about the Z axis instead of 
the Y one. 

As a example, consider what Qubiter gives for the Discrete Fourier Transform 
matrix U. This matrix has entries 



where a, b G BooI Nb . N s and d(-) were defined in Section 2(a). When N B = 2, 
Qubiter gives 33 operations (see Fig. 8). After doing the trivial optimization of 
removing all factors for which the rotation angle is zero, the 33 operations in Fig.8 
reduce to 25 operations in Fig.9. In the future, we plan to introduce into Qubiter 
many more optimizations and some quantum error correction. Much work remains 
to be done. 




1 



exp 



i2nd(a)d(b) 



(5.1) 
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Fig.8 Output of Qubiter with zero angle optimization turned OFF. (input: 2-bit 

Discrete Fourier Transform matrix). 
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ROTZ 292.500000 

ROTZ 1 45.0000000 

CNOT 1 T 

ROTZ 315.000000 

CNOT 1 T 

ROTY 315.000000 

ROTZ 337.500000 

CNOT 1 T 

ROTZ 45.0000000 

CNOT 1 T 

ROTY 1 315.000000 

CNOT T 1 

ROTY 1 22.5000000 

CNOT T 1 

ROTZ 247.500000 

ROTZ 1 45.0000000 

CNOT 1 T 

ROTZ 315.000000 

CNOT 1 T 

ROTY 315.000000 

ROTZ 22.5000000 

CNOT 1 T 

ROTZ 315.000000 

CNOT 1 T 

PHAS 157.500000 

Fig. 9 Output of Qubiter with zero angle optimization turned ON. (input: 2-bit 

Discrete Fourier Transform matrix). 
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FIGURE CAPTIONS: 

FlG.l A binary tree. Each node (3 has a single parent. If the parent is to /?'s right 
(ditto, left), then (3 contains the names of the matrices produced by applying 
the CS Decomposition Theorem to the L matrices (ditto, R matrices) of /?'s 
parent. 



Fig. 


.2 


Circuit symbols for the 4 different types of c-nots. 


Fig. 


3 


Four equivalent circuit diagrams for Exchanger. 


Fig. 


.4 


Four equivalent circuit diagrams for Twin-to-twin-er. 


Fig. 


.5 


Circuit symbol for Exchanger. 


Fig. 


.6 


Circuit diagram for Eq.(3.12). 


Fig. 


.7 


Circuit diagram for Eq.(3.14). 


Fig. 


.8 


Output of Qubiter with zero angle optimization turned OFF. (input 



Discrete Fourier Transform matrix). 

FlG. 9 Output of Qubiter with zero angle optimization turned ON. (input: 2-bit 
Discrete Fourier Transform matrix). 
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